Comparative transcriptome analysis reveals the biocontrol mechanism of Bacillus velezensis E68 against Fusarium graminearum DAOMC 180378, the causal agent of Fusarium head blight

Fusarium graminearum is the causal agent of Fusarium Head Blight, a serious disease affecting grain crops worldwide. Biological control involves the use of microorganisms to combat plant pathogens such as F. graminearum. Strains of Bacillus velezensis are common biological control candidates for use against F. graminearum and other plant pathogens, as they can secrete antifungal secondary metabolites. Here we study the interaction between B. velezensis E68 and F. graminearum DAOMC 180378 by employing a dual RNA-seq approach to assess the transcriptional changes in both organisms. In dual culture, B. velezensis up-regulated genes related to sporulation and phosphate stress and down-regulated genes related to secondary metabolism, biofilm formation and the tricarboxylic acid cycle. F. graminearum up-regulated genes encoding for killer protein 4-like proteins and genes relating to heavy metal tolerance, and down-regulated genes relating to trichothecene biosynthesis and phenol metabolism. This study provides insight into the molecular mechanisms involved in the interaction between a biocontrol bacterium and a phytopathogenic fungus.


Introduction
Fusarium graminearum is the causal agent of Fusarium Head Blight (FHB) in wheat and barley [1]. It is considered a major global threat that impacted great economic losses on the cereal industry because of reduced grain yield and grain quality. FHB can also cause contamination of crops and grains with diverse mycotoxins, including deoxynivalenol (DON) and zearalenone, which are harmful for humans and animals [2,3]. FHB has severe impacts on grain yield, with reported losses of up to 80% [4]. It has been estimated that in 2015, FHB costed the US wheat and barley industries $1.469 billion in value of yield forgone [5]. Recent data for Canada is not available, though it has been estimated that $520 million were lost by wheat graminearum in vitro and limit FHB disease progression in planta when applied as biocontrol agents [22][23][24]. The bacterial strains from these studies and similar ones are typically endophytes isolated directly from plants, ensuring their ability to colonize the hosts. However, it is important to assess the abilities of bacterial strains from other environments as they may have differing effects as biocontrol candidates as compared to endophytes. Bacillus velezensis strain E68 was previously isolated from crude oil samples from the Albertan oil sands [25,26]. Oil-dwelling bacteria must be adapted to live in the harsh conditions: high temperature, high pressure and low carbon [27,28]. These bacteria have adapted to the particular organic carbon sources present in these environments, including hydrocarbons. Biosurfactants secreted by these bacteria are used to degrade hydrocarbons for use as a carbon source [29]. Much research has gone into studying the potential of these biosurfactants for bioremediation of hydrocarbon pollutants [30]. However, oil-dwelling bacteria are unexplored when it comes to biocontrol. In preliminary experiments, B. velezensis E68 was found to inhibit the growth of various fungal phytopathogens such as Rhizoctonia solani and Fusarium graminearum (S1 Fig). Also, the cell-free supernatant of strain E68 exhibited biosurfactant activity using drop collapse and agar spreading assays [26]. These results suggest strain E68 has a wide inhibition spectrum against plant fungal pathogens.
The main task of the present research is to dissect the undergoing changes in the transcriptome profiles of both partners following dual culture of B. velezensis E68 and F. graminearum DAOMC 180378 using RNA-seq. This could provide valuable insights into the main factors that determine its outcome. Using the RNA-seq method, we observed the transcriptional activity of both organisms in response to the presence of one another.

Biological material and culture conditions
Cultures of Fusarium graminearum strain DAOMC 180378 were provided by the Canadian Collection of Fungal Cultures (DAOMC; Agriculture and Agri-Food Canada). Bacillus velezensis strain E68 was previously isolated from an oil well in Alberta [25,26]. Microbial strains were retrieved from 20% glycerol stock stored at -80˚C. F. graminearum DAOMC 180378 was grown on potato dextrose agar (PDA; Difco Laboratories) at 25˚C in the dark. B. velezensis E68 was grown on Luria-Bertani agar (LBA) at 37˚C.

Inhibitory effects of Bacillus velezensis E68 on Fusarium graminearum DAOMC 180378 mycelial growth
In order to test the inhibitory effect of Bacillus velezensis E68 on Fusarium graminearum DAOMC 180378 mycelial growth, the two organisms were grown in dual culture. Sterilized circular cellophane membranes (UCB, North Augusta, USA) were overlaid on 9 cm diameter petri plates containing 20 mL of PDA. Two half-moon sections were cut out from opposite sides of the membrane up to 2 cm from its center in order for B. velezensis to grow directly on the culture media (Fig 1). An agar plug (5 mm diameter) was taken from 7-day old cultures of actively growing F. graminearum, and were placed at the center of each plate overlaid with cellophane. A 10 μL aliquot of Bacillus velezensis (OD 600 of 1.0) from an overnight LB broth was spotted on 2 opposite edges of the plate where the cellophane membrane was cut out, 2.5 cm from the center of the plate. A 0.5 cm agar plug of F. graminearum placed at the center of a PDA plate overlaid with a cellophane membrane served as control (Fig 1). Three 10 μL aliquots of Bacillus velezensis (OD 600 of 1.0) spotted apart on a PDA plate served as control (Fig 1). Dual culture and control plates were incubated at 24˚C for 3 days before RNA extraction.
Dual culture and control plates were set up for microscopy as described above. A strip of cellophane membrane (0.5 cm x 2 cm) containing edges of the fungal hyphae was carefully excised from both a dual culture and a control treatment plate. The strips were placed on a microscopy slide, stained with bromophenol blue and mounted for microscopy on an EVOS XL Core Imaging System (Invitrogen).

Sample preparation and RNA extraction
Bacterial colonies were collected with a sterile loop by removing the half of each colony that was facing the center of the plate (Fig 1). In total, 6 halves of bacterial colonies were pooled into one RNA sample (i.e. 3 dual-culture plates pooled into 1 sample and 1 control plate pooled into 1 sample). The bacterial cells was placed into 500 μL of sterile water per sample for immediate extraction. Strips of fungal hyphae (5 x 50 mm) closest to the bacterial spots were collected with a sterile spatula (Fig 1). Six hyphal strips were pooled into one RNA sample. For the control treatment, an equivalent amount of hyphae was taken from the edges of a fungal colony and used as one sample. The fungal material was immediately ground with a mortar and pestle in liquid nitrogen and stored at -80˚C before RNA extraction.

PLOS ONE
High quality total RNA from B. velezensis and F. graminearum was extracted using the Nucleospin RNA kit (Takara). RNA extraction from B. velezensis followed manufacturer's instructions for Gram-positive bacteria, with modifications: the cell pellet was resuspended in 100 μL TE buffer (10 mM Tris, 1 mM EDTA, pH 8.0) with 3 mg/mL lysozyme, then incubated for 15 minutes at 37˚C. RNA extraction from F. graminearum followed the instructions for cultured cells and tissues. For both organisms, 600 μL of buffer RA1, 6 μL of β-mercaptoethanol and 600 μL of 70% ethanol were used instead of 350 μL, 3.5 μL and 350 μL respectively. In all cases, elution was done by eluting twice in 60 μL of RNase-free water. RNA samples integrity and concentration were assessed by 1% bleach agarose gel and Nanodrop ND-1000. Samples were stored at -80˚C for downstream applications.

RNA sequencing
Five RNA biological replicates for each treatment (B. velezensis control and dual-culture, F. graminearum control and dual-culture) were sent for library preparation and sequencing at IDSeq Inc, (Sacramento, California). The B. velezensis libraries were prepared with a TruSeq RNA library prep kit, with a Ribo-Zero Plus rRNA depletion kit. The F. graminearum libraries were prepared with the NEBNext Ultra II RNA library prep kit with poly-A selection. Samples were sequenced on an Illumina HiSeq 4000.

RNA-seq analysis
Fig 2 depicts the pipeline followed from sample preparation to RNA-seq data analysis. Sequencing reads were checked for quality and adaptor contamination before and after trimming with FastQC 0.11.9 [31]. Low quality bases and adaptor sequences were trimmed from the raw sequencing reads with Trim Galore 0.6.5 with stringency set to 3 [32]. The genome assembly for F. graminearum DAOMC 180378 is available on GenBank, however no reference annotation is available. The genome annotation version 5 [33] from Fusarium graminearum PH-1 was transferred to the genome assembly of DAOMC 180378 using Liftoff 1.6.1 [34]. The resulting annotation was converted from GFF3 format to GTF format with gffread 0.12.3 [35]. Trimmed reads were mapped to the reference genomes for F. graminearum DAOMC180378 and Bacillus velezensis E68 using STAR 2.7.5a on default settings [25,36,37]. The GeneCounts function from STAR was used to obtain read counts. Differential gene expression was assessed using edgeR 3.32.1 [38]. Genes were filtered by CPM expression using edgeR default settings. The edgeR GLM framework was used. Genes with a log 2 fold change under 1.5 were deemed non-significant using the glmTreat function. Genes with a p-value under 0.05 after correction for multiple testing were deemed significant and differentially expressed. edgeR was also used to perform rotation gene set testing (ROAST) [39]. Genes were functionally annotated with gene ontology data by using Trinotate 3.2.2 to facilitate BLASTx searches of the gene coding sequences against the UniProtKB protein database [40]. Gene ontology enrichment analysis was carried out with GOSeq 1.44.0 to detect enriched gene ontology terms among up and down-regulated genes in both organisms, using a false discovery rate < 0.05 [41]. Biosynthetic gene clusters in B. velezensis E68 were identified with antiSMASH in previous work [25].

Validation of RNA-seq by qRT-PCR
In order to confirm the differential expression results obtained by RNA-seq, 8 genes for each organism were selected for validation by qRT-PCR (Table 1). Three RNA samples from each condition were used to validate the RNA-seq. The RNA samples were transcribed into cDNA using the Agilent AffinityScript qPCR cDNA Synthesis Kit. Specific primers for the selected genes were designed using Geneious 8.1.9 and the OligoAnalyzer tool from IDT. Primers are Cultures of B. velezensis E68 and F. graminearum DAOMC 180378 were grown in single and dual culture. The striped areas were collected with a sterile spatula before RNA extraction. Total RNA was sent for library preparation and RNA sequencing. Bacterial RNA was depleted for ribosomal RNA, while fungal RNA was enriched with poly-A selection. Samples were sequenced on an Illumina HiSeq using a paired-end 150 bp protocol. Raw sequencing reads were preprocessed by removing adaptors and poor quality bases and reads. Filtered reads were aligned to their respective genomes and counted for each gene. Gene counts were normalized by the TMM method, low expression genes were removed from analysis and differential expression was calculated. In order to functionally annotate the genes, each gene sequence was used in a BLASTx search to compare to the UniprotKB protein database, annotating each gene with a description and gene ontology terms. The resulting GO terms were used in GO enrichment analysis. Known biosynthetic gene clusters were used in rotation gene set testing to determine their overall regulation pattern. https://doi.org/10.1371/journal.pone.0277983.g002

PLOS ONE
listed in Table 1. Specificity was assessed with PCR and agarose gel. recA (GYA98_RS07520) and secA (GYA98_RS15380) were used as housekeeping genes for B. velezensis [42]. GzUBH, ubiquitin thiolesterase (FGRAMPH1_01G03049) and TUB2, tubulin beta chain (FGRAMPH1_01G26865) were used as housekeeping genes for F. graminearum [43,44]. Dilutions of pooled cDNA were used to create standard curves for each primer to calculate reaction efficiency and determine the optimal cDNA concentration for final usage. The geometric mean of the housekeeping genes was used as a normalization factor. The reaction efficiency and normalization factor were used to adjust the ΔC q values of the genes of interest as described in Taylor et al. [45]. Log transformed normalized expression values were compared between single and dual culture using a t-test. Three technical replicates were run per sample.

Dual RNA-seq analysis of B. velezensis and F. graminearum
Cultures of Bacillus velezensis E68 and Fusarium graminearum DAOMC 180378 were grown together on potato dextrose agar as described in the Methods. Following 3 days of incubation, a clear inhibition of fungal growth was observed ( Fig 3B). The morphology of F. graminearum appeared healthy in monoculture ( Fig 3C). When grown in dual culture with B. velezensis, the morphology of F. graminearum appeared collapsed, with hyphal tips disrupted and developing circular spheres (Fig 3D and 3E). To evaluate the changes in total gene expression in both organisms, RNA was separately extracted from both organisms growing alone in single and together in dual culture. RNA from both organisms from both conditions was sequenced.
After adaptor trimming and filtering the sequencing reads for quality, between 2.9 and 10.4 million sequencing reads were obtained for B. velezensis samples and between 16.4 and 24.5 million reads were obtained for F. graminearum samples ( Table 2). The remaining reads were aligned to the reference genomes and reads aligning to each gene were quantified. The R package edgeR was used to filter out lowly expressed genes, normalize the libraries and calculate differential expression for each gene between single and dual culture. A total of 168 genes from B. velezensis and 4207 genes from F. graminearum were filtered out due to low expression. After calculation of differential expression, genes with a log 2 fold change under 1.5 between conditions were not considered significant, as well as genes with a false discovery rate above 0.05. With these requirements, 1155 genes (32.4% of all expressed genes) were differentially expressed in B. velezensis between single and dual culture. In F. graminearum, 1503 genes (15.3% of all expressed genes) were differentially expressed (Table 3). Of the differentially expressed genes (DEGs), 47.2% were down-regulated in B. velezensis, whereas 66.1% of DEGs were down-regulated in F. graminearum. Genes were functionally annotated by BLASTx searches to the UniprotKB database, which also annotated many genes with gene ontology (GO) terms. In B. velezensis, 2394 genes were functionally annotated with GO terms (67.2% of all expressed genes), while 5464 genes were functionally annotated with GO terms in F. graminearum (55.6% of all expressed genes).

Transcriptomic changes of B. velezensis in dual culture with F. graminearum
Following gene ontology enrichment analysis, 16 GO terms were enriched among significantly up-regulated genes in B. velezensis, whereas 7 GO terms were enriched among significantly down-regulated genes (Table 4). Genes associated to sporulation and related GO terms such as the spore wall and spore germination were up-regulated. Many genes were annotated as spore coat proteins as well as many different genes from sporulation stages II, III, IV and V. These included spoIIE, spoIIGA, spoIIR, and sigE from stage II, spoIID, spoIIP, spoIIQ, spoIIIAH, FisB, sigG and sigK from stage III, spoIVB, ctpB, spoIVA, sleB and cwlJ from stages IV and V (S2 Table). The master sporulation regulator gene spo0A which initiates sporulation (locus tag GYA98_RS10380) was not significantly differentially regulated (S1 Table). The phosphotransferase spo0B (GYA98_RS11640) was also not regulated, but the phosphotransferase spo0F (GYA98_RS16270) was down-regulated. The GO term for oxidoreductase activity, acting on CH or CH2 groups, was enriched in up-regulated genes. Genes associated with this term included 2 subunits of ribonucleoside-diphosphate reductase: nrdE (GYA98_RS07750) and nrdF (GYA98_RS07755) (S1 Table). This enzyme contributes to DNA synthesis by converting ribonucleotides to deoxyribonucleotides. The 4 other up-regulated genes associated to this term were subunits of xanthine dehydrogenase (GYA98_RS17550, 17555, 17560, 17570) (S1 Table). Another subunit GYA98_RS17565 was also up-regulated, though not annotated with the GO term. Xanthine dehydrogenase participates in purine degradation by degrading xanthine into uric acid. Two terms for teichuronic acid biosynthesis were enriched in up-regulated

PLOS ONE
genes (Table 4). All genes associated to these terms were part of the tua operon. These genes synthesize teichuronic acid, which replaces teichoic acid in the cell wall in phosphate starvation conditions. In addition, all 5 genes annotated with the ATPase-coupled phosphate ion transmembrane transporter activity GO term were up-regulated. These genes are all part of  Table). Both the tua and PstSABC operons are controlled by the Pho regulon. Another member of this regulon, glpQ (GYA98_RS00270), was also up-regulated in dual culture [47] (S2 Table). Secondary metabolites in B. velezensis are produced through enzyme complexes encoded by biosynthetic gene clusters. The transcription of these gene clusters was assessed by treating each cluster as a gene set and using rotation gene set testing (ROAST). The genome of B. velezensis E68 contains biosynthetic gene clusters which may synthesize 9 known secondary metabolites. Rotation gene set testing showed that all gene clusters were significantly downregulated, with the exception of surfactin, which was not significant (Table 5). This was also reflected in gene ontology, as GO terms related to antibiotic biosynthesis and metabolism were enriched among down-regulated genes in dual culture (Table 4). Despite being down- regulated, genes responsible for synthesis of secondary metabolites remained fairly highly expressed in dual culture. The median transcripts per million (TPM) of all genes across all dual culture B. velezensis samples was 73.207. In comparison, the expression of genes annotated as core biosynthetic genes by antiSMASH remained high in down-regulated clusters encoding secondary metabolites (Table 6). GO terms related to the tricarboxylic acid (TCA) cycle were also enriched among downregulated genes. Significantly down-regulated genes associated to these terms included the   Table). Malate dehydrogenase (GYA98_RS12230), which was not annotated with the TCA cycle GO terms, was also significantly down-regulated. The GO term for polysaccharide biosynthetic process was also enriched among the down-regulated genes. Of the significantly down-regulated genes associated to this term, 14 of the 15 genes of the eps operon (epsB-O, GYA98_RS14870-14935) were included (S2 Table). The epsA homolog in B. velezensis E68 (GYA98_RS14940) was not significantly differentially expressed. The eps operon encodes enzymes responsible for synthesis of exopolysaccharide, a component of the extracellular biofilm matrix [48,49]. In addition, tasA (GYA98_RS10575), the major protein component of the extracellular matrix, was also downregulated in dual culture, though it wasn't annotated with GO terms.

Transcriptomic changes of F. graminearum in dual culture with B. velezensis
In F. graminearum, no GO terms were significantly enriched among up-regulated genes in dual culture. Despite not being significantly enriched, 5 GO terms had all 3 of their associated genes be up-regulated in F. graminearum (Table 4). FGRAMPH1_01G11895 and FGRAMPH1_01G06603, both annotated as heavy metal tolerance proteins which may transport phytochelatin, and FGRAMPH1_01G26673, an ATP-binding cassette transporter, were all up-regulated in dual culture conditions and were the only genes associated with the 5 GO terms for phytochelatin transport/import and glutathione import into the vacuole (S3 Table). The relative lack of total genes associated to the terms prevented these terms from being significantly enriched after correction for multiple testing. The genome of F. graminearum DAOMC 180378 was annotated with 4 genes described as killer protein 4-like. All 4 genes were up-regulated in dual culture (S4 Table).
In the 16 GO terms enriched among down-regulated genes in F. graminearum, many were related to the broad categories of the catabolic process term and the oxidoreductase activity term. More specifically, many of the terms related to the metabolism of aromatic compounds, notably the terms for phenol containing, benzene-containing and catechol-containing compound metabolic processes, as well as the term for 3,4-dihydroxybenzoate (also known as protocatechuate) metabolic process (Table 4). Significantly down-regulated genes associated to these compounds included genes upstream of protocatechuate/catechol: 2 genes annotated as quinate dehydrogenases (FGRAMPH1_01G13935 and FGRAMPH1_01G18555), one gene annotated either as a quinate or shikimate dehydrogenase (FGRAMPH1_01G13945),  Table). In addition, genes degrading catechol and protocatechuate were down regulated: a catechol 1,2-dioxygenase (FGRAMPH1_01G13423) and another gene annotated as a catechol 1,2-dioxygenase, but a BLAST search against the DAOMC 180378 genome showed that protocatechuate-3,4-dioxygenase from Aspergillus nidulans was homologous to this gene (FGRAMPH1_01G12667) [50]. Also part of the significantly down-regulated genes associated to benzene compounds were 3 genes from the homogentisate pathway, which utilizes phenylacetate and phenylalanine as carbon sources. Homogentisate 1,2-dioxygenase (FGRAMPH1_01G05639), 3-hydroxyphenylacetate 6-hydroxylase (FGRAMPH1_01G22299) and fumarylacetoacetase (FGRAMPH1_01G11459) were all down-regulated (S3 Table). In addition, 3 genes from the kynurenine pathway were down-regulated: 2 kynureninases (FGRAMPH1_01G16409 and FGRAMPH1_01G14561), as well as an indoleamine 2,3-dioxygenase (FGRAMPH1_01G16407) (S3 Table). This pathway degrades tryptophan into NAD+. In dual culture, F. graminearum hyphae was distorted and swelling occurred, based on microscopy (Fig 3). However, the differential gene expression data from RNA-seq did not show concerted differential regulation of genes relating to this phenotype. This included genes related to ergosterol synthesis, glucan or chitin synthases and genes related to autophagy. Similarly, F. graminearum also did not show large scale differential regulation of genes relating to oxidative stress, specifically genes encoding for superoxide dismutases and catalases (S4 Table). Genes annotated as glutathione S-transferases were also mostly not differentially expressed, though some were significantly down-regulated in dual culture (S4 Table).
F. graminearum DAOMC 180378 is known to produce deoxynivalenol, a trichothecene mycotoxin [37]. To assess the effect of dual culture with B. velezensis on trichothecene mycotoxin expression, genes from the core trichothecene gene cluster, as well as other genes from the trichothecene pathway were evaluated for differential expression (Table 7) [51]. Biosynthetic enzymes TRI3, TRI8 and TRI11 did not show significant differential expression. The biosynthetic enzymes TRI1, TRI5, TRI13 and TRI101 were significantly down-regulated, while TRI4 was down-regulated, only approaching significance. The accessory proteins TRI6, TRI9,

PLOS ONE
TRI10, TRI11 and TRI14 did not show significant differential expression. The transcription factor TRI15 showed significant up-regulation in dual culture. TRI15 is likely a negative regulator of trichothecene expression [52].

Validation of RNA-seq via quantitative PCR
In order to validate the results of RNA-seq, 8 genes were selected per organism to evaluate their differential expression through qRT-PCR (Table 1). Three RNA samples per condition were used for analysis (Bv1, Bv2, Bv5, DB1, DB2, DB5, Fg1, Fg2, Fg5, DF1, DF2, DF5). The housekeeping genes were selected based on use in previous studies as well as lack of differential expression in the RNA-seq dataset. In B. velezensis E68, the 2 selected non-differentially expressed genes (hisG and sigM) saw identical results in qRT-PCR (Fig 4A). The gene hssR showed significant up-regulation in both RNA-seq and qRT-PCR. Four significantly differentially expressed genes (phoD, sigW, cotE, citZ) had similar trends between RNA-seq and qRT-PCR, though the qRT-PCR fold-changes were not significant after t-test. One gene, catE showed opposite trends in qPCR, though the result was not significant. In F. graminearum DAOMC 180378, the 2 selected non-differentially expressed genes RVB2 and RAD5 showed significant differential expression using qRT-PCR (Fig 4B). The 3 genes YBT1, AES1 and DMLA showed significant differential expression in the same direction in both RNA-seq and qRT-PCR. The genes ZRA1 and GRA1 showed significant regulation in opposite directions in RNA-seq and qRT-PCR. The ALDH gene also showed opposite regulation by RNA-seq and qRT-PCR, though the qRT-PCR result was not significant.

Discussion
A dual RNA-seq approach was developed to compare the transcriptional profiles of two organisms, Bacillus velezensis E68 and Fusarium graminearum DAOMC 180378 in single and dual culture. Previously, dual RNA-seq was employed to investigate the interaction between the rhizobacteria Lysobacter capsici and the oomycete Phytophthora infestans [53]. In addition, a dual microarray approach was used to assess the transcriptional changes from both sides of a Bacillus subtilis and Aspergillus niger interaction [54]. However, to the best of our knowledge, this is the first time a dual RNA-seq approach has been used to study the transcriptional changes between a biocontrol candidate bacteria and a filamentous fungus. Strains of B. velezensis and the closely related B. subtilis were previously shown to inhibit growth of F. graminearum on solid media, demonstrating the potential of these strains as biological control candidates [23,24,55]. In our study, B. velezensis E68 showed very similar inhibition of F. graminearum DAOMC 180378. Microscopic observations of F. graminearum exposed to B. velezensis showed distorted hyphal tips and circular swellings, most likely due to secretion of antifungal lipopeptides by B. velezensis. This phenotype was also observed in other studies exposing F. graminearum to both B. velezensis cultures and also purified lipopeptide extracts from B. velezensis and B. subtilis [23,[56][57][58][59]. Lipopeptides secreted from B. velezensis and similar species are known to cause damage to the plasma membrane and cell wall of various fungal phytopathogens, including F. graminearum [23,58,60,61]. When confronted with F. graminearum, B. velezensis showed clear up-regulation of many genes related to sporulation. Strains of Bacillus tend to sporulate under conditions of limited nutrient availability, though the specific signals by which this is determined are unknown [62]. In our study, genes related to stages 0 and I of sporulation did not show significant differential expression in B. velezensis. This includes the master regulator of sporulation spo0A, though since its activity is reliant on its phosphorylation state, its transcriptional status may not be of importance [63]. However, multiple genes from stages II through V were up-regulated. SigF, the first sigma factor related to sporulation is produced prior to stage II. The gene product of spoIIE is responsible for activating this sigma factor in the forespore [64]. Consistent with this idea, in our study spoIIE was up-regulated in B. velezensis E68 in dual culture. In addition, spoIIGA and spoIIR, the activators of sigE, as well as sigE itself, were all up-regulated in dual culture, triggering stage III: engulfment [65]. Genes required for engulfment and membrane fission were up-regulated in dual culture: spoIID, spoIIP, spoIIQ, spoIIIAH and FisB [66,67].

PLOS ONE
The forespore-specific sigma factor sigG and the mother cell specific sigma factor sigK are both activated at the end of engulfment, these are both up-regulated in dual culture. In addition the sigK activators spoIVB and ctpB were up-regulated [68][69][70]. Many spore cortex coat proteins related to stages IV and V were up-regulated, including spoIVA, which forms the base layer of the spore coat [68]. Furthermore, spore germination genes were also up-regulated in dual culture. sleB and cwlJ are enzymes which lyse peptidoglycan during germination and were both up-regulated. However, they lie dormant in the spore in a mature form [71]. Taken together, these results suggest that in the presence of F. graminearum, B. velezensis cells sporulate, at stages II through V of sporulation on day 3 during our dual culture assay. Previous research showed B. subtilis up-regulated late sporulation genes when grown on Aspergillus niger hyphae, though no significant difference in spore content was found [54].
In dual culture, B. velezensis E68 down-regulated key components of the extracellular biofilm matrix: the eps operon and the tasA gene [72]. The decision to begin biofilm formation and the decision to sporulate are both governed by the Spo0A protein [73]. Spo0A is a transcriptional regulator whose activity depends on its phosphorylation state. High levels of phosphorylated Spo0A promote entry into sporulation through repression of sinI expression [63].
Research also shows that sporulation may require biofilm formation to occur [73]. As a result, the down-regulation of biofilm components in B. velezensis coincides with the up-regulation of sporulation genes.
The genes encoding for the tricarboxylic acid cycle enzymes were down-regulated in B. velezensis in dual culture. This is consistent with previous research studying the transcriptomics of the life cycle of Bacillus pumilus [74]. B. pumilus showed decreased expression of TCA cycle genes in the stationary growth phase as compared to the exponential growth phase, in addition to an increase in expression of sporulation-related genes, much the same as B. velezensis in our experiment. Furthermore, TCA cycle genes were down-regulated and sporulation genes were up-regulated in B. subtilis which were attached to Aspergillus niger hyphae as compared to planktonic cells [54].
Secondary metabolism in B. velezensis E68 was evidently altered, with all known biosynthetic gene clusters (with the exception of surfactin) showing reduced expression in dual culture with F. graminearum. Despite the down-regulation, these metabolites remained expressed in dual culture and likely secreted into the media, as reflected by the changes in fungal hyphae appearance in dual culture (Fig 3). This down-regulation is unexpected, as multiple studies of lipopeptide-secreting Bacillus report up-regulation of lipopeptide related genes or increased lipopeptide concentrations in confrontation with fungi [75][76][77][78][79]. However, surfactin was found to be down-regulated in B. subtilis cells attached to hyphae of A. niger [54]. Our experiment was performed at a single time-point: 3 days after establishment of interaction. It may be possible that different regulation patterns could be expressed for secondary metabolism genes in earlier or later stages of the interaction. Using proteomic methods for the detection of lipopeptides and other secondary metabolites could also show variations from the gene expression data obtained here.
In dual culture, B. velezensis also showed signs of phosphate stress. In times of phosphate starvation, Bacillus strains induce genes of the Pho regulon, controlled by the PhoP sensor and the PhoR regulator, though these were not differentially expressed in our study [47]. When confronted with F. graminearum, B. velezensis up-regulated members of the Pho regulon, such as the tua operon. This operon synthesizes teichuronic acid, a component of the cell wall. In phosphate starvation conditions, teichuronic acid replaces teichoic acid in the cell wall [80]. Phosphate can then be collected from the teichoic acid through the action of phoD, which was also up-regulated in dual culture [81]. The PstSABC operon is another part of the Pho regulon which was up-regulated in dual culture. The PstSABC complex is a phosphate importer which may also play a role in sensing extracellular phosphate concentrations through interaction with PhoR [46,82]. Furthermore, B. velezensis up-regulated glpQ, also part of the Pho regulon. This gene encodes a glycerophosphoryl diester phosphodiesterase, which hydrolyzes deacetylated phospholipids [83].
Taken together, these results indicate that in dual culture with F. graminearum, B. velezensis E68 entered sporulation, as shown by the up-regulation of sporulation and biofilm formation genes. This could be due to different stresses, including phosphate starvation. Consistent with its entry into sporulation, B. velezensis down-regulated genes related to carbohydrate metabolism and secondary metabolite biosynthesis.
In our study, there were no GO terms enriched among up-regulated genes in F. graminearum. This can be partially attributed to the low number of expressed genes annotated with GO terms (55.6% of expressed genes). Additionally, only 510 genes were significantly up-regulated, as compared to 993 down-regulated genes. As a result, global up-regulation of particular molecular functions or biological processes was not found. When the fungal phytopathogen Sclerotinia sclerotiorum was confronted to B. velezensis fermentation broth, the fungus was found to induce genes related to cell wall and ergosterol synthesis, fatty acid synthesis, antioxidants and autophagy [84]. In this study, genes related to cell wall synthesis via either glucan synthase or chitin synthase were mostly not differentially expressed, with chitin synthase 6 (FGRAMPH1_01G22563) being the only significantly up-regulated chitin synthase and chitin synthase 5 (FGRAMPH1_01G10489) being down-regulated. Distinct up-regulation of ergosterol pathway genes as defined through annotation or previous evidence was also not found. It may be that fungal cells most affected by the bacterial lipopeptides were lysed and their RNA could not be sequenced. Genes related to autophagy were also not differentially expressed in dual culture [85]. In addition to damage to the cell wall and membranes, Bacillus lipopeptides can also induce accumulation of reactive oxygen species [58,86]. However, our data signalled that F. graminearum did not show concerted up-regulation of genes relating to oxidative stress. Five genes previously identified as superoxide dismutases did not show differential expression [87]. From a list of 7 previously identified catalases in F. graminearum, 6 showed no differential expression, while 1 was significantly down-regulated [88]. One catalase which was separately identified using its annotation was up-regulated in dual culture (FGRAMPH1_01G05343). In addition, there were no genes annotated as glutathione S-transferase that were up-regulated in dual culture, though several were significantly down-regulated. F. graminearum did show up-regulation of 3 genes related to heavy metal tolerance through phytochelatin or glutathione transport into the vacuole. Phytochelatin is a polymer of glutathione which is used to sequester cadmium. Upon transport into the vacuole via HMT1, an ABC transporter, the phytochelatin-cadmium complex can integrate sulfide to achieve higher cadmium binding and stability. Phytochelatin may also confer tolerance to copper, mercury, silver and arsenate ions [89].
While not associated to particular GO terms, all 4 genes previously identified as killer protein 4-like proteins were up-regulated in dual culture. These had previously been shown to be induced by F. graminearum in stressful conditions such as high salt concentrations or osmotic pressures, and it is believed they play a role in virulence, though it is uncertain which specific function they perform in the stress response [90].
F. graminearum showed a clear down-regulation of genes associated to metabolism of phenolic compounds in dual culture with B. velezensis. In nature, F. graminearum overwinters on infested crop residues, adopting a saprotrophic lifestyle [91]. Fungi living on plant litter must decompose the lignin-related aromatic compounds to use as a carbon source [92]. Fungi employ various pathways to degrade and recycle aromatic compounds, usually through one of 5 intermediates: protocatechuate, catechol, homogentisate, gentisate or hydroquinone [50].
Potato, a major component of the culture media, is known to contain various phenolic compounds such as phenolic acids and lignin [93,94]. F. graminearum showed down-regulation of genes in the protocatechuate, catechol and homogentisate pathways. It also down-regulated genes of the kynurenine pathway, responsible for degrading the aromatic amino acid tryptophan into NAD+. It is possible that in dual culture, B. velezensis is metabolizing some of the aromatics in the media, lessening the requirements for aromatic metabolism in the fungus. Phenolic compounds such as 2,3-dihydroxybenzoate are required for the synthesis of bacillibactin, a secondary metabolite produced by B. velezensis [95]. This notion merits further investigation.
In dual culture, F. graminearum showed significant down-regulation of 3 biosynthetic enzymes: TRI1, TRI5 and TRI13, with TRI4 approaching significance. TRI1, TRI4 and TRI5 are all key components of the trichothecene synthesis pathway [51]. TRI4 and TRI5 carry out the first 5 steps of the biosynthetic pathway, while TRI1 performs a later step. TRI101 was significantly up-regulated in dual culture, however, which is unexpected given its position between TRI4 and TRI1 in the biosynthetic pathway. TRI13 has no active functionality in deoxynivalenol-producing chemotypes of F. graminearum, which include strain DAOMC 180378 [51]. This is reflected by its low expression relative to the other genes of the TRI locus. TRI15 is a negative regulator of trichothecene biosynthesis which was up-regulated in dual culture [52]. Its specific gene targets are unknown, though exogenously added trichothecenes induced its expression in F. sporotrichoides, a trichothecene producer. Studies of trichothecene gene expression have not analyzed TRI15 expression, though we speculate its up-regulation in dual culture caused down-regulation of trichothecene biosynthetic genes. Down-regulation of trichothecene biosynthetic genes or reduction in trichothecene production through dual culture with Bacillus has been observed before [22][23][24]. The mechanism through which this gene inhibition occurs remains to be elucidated.
In summary, this work investigated the transcriptomic response of both a biocontrol candidate bacteria and its fungal phytopathogenic target in a dual culture assay. B. velezensis E68 was found to up-regulate genes related to sporulation and phosphate stress response and down-regulate genes related to secondary metabolism, biofilm formation and the TCA cycle. F. graminearum DAOMC 180378 saw up-regulation of heavy metal tolerance genes and killer protein-4 like proteins. It down-regulated trichothecene biosynthesis and phenol metabolism.
Supporting information S1